home *** CD-ROM | disk | FTP | other *** search
/ SGI Developer Toolbox 6.1 / SGI Developer Toolbox 6.1 - Disc 4.iso / public / SciAn / src / ScianVisSticks.c < prev    next >
C/C++ Source or Header  |  1994-08-01  |  33KB  |  1,066 lines
  1. /*
  2.   ScianVisSticks.c
  3.   Tzong-Yow Hwu
  4.   Modified from ScianVisBalls.c
  5.   Sticks visualization object in SciAn  
  6. */
  7. /* 
  8.    Needs to do:
  9.    1.  the A-corbon chain dealing way.
  10.    2.  joint.
  11. */
  12.  
  13. #include "Scian.h"
  14. #include "ScianTypes.h"
  15. #include "ScianArrays.h"
  16. #include "ScianWindows.h"
  17. #include "ScianTextBoxes.h"
  18. #include "ScianButtons.h"
  19. #include "ScianTitleBoxes.h"
  20. #include "ScianObjWindows.h"
  21. #include "ScianIcons.h"
  22. #include "ScianColors.h"
  23. #include "ScianControls.h"
  24. #include "ScianLists.h"
  25. #include "ScianSliders.h"
  26. #include "ScianIDs.h"
  27. #include "ScianDatasets.h"
  28. #include "ScianPictures.h"
  29. #include "ScianDialogs.h"
  30. #include "ScianErrors.h"
  31. #include "ScianComplexControls.h"
  32. #include "ScianMethods.h"
  33. #include "ScianStyle.h"
  34. #include "ScianVisObjects.h"
  35. #include "ScianVisSticks.h"
  36. #include "ScianDraw.h"
  37. #include "ScianTemplates.h"
  38. #include "ScianTemplateHelper.h"
  39. #include "ScianSymbols.h"
  40.  
  41. ObjPtr visSticks;              /*Class for Sticks objects*/
  42.  
  43. /* Stick Styles */
  44. #define BAS_LINES     0
  45. #define BAS_GIRDERS   1
  46. #define BAS_CYLINDERS 2
  47.  
  48. /* color shading definition */
  49. #define BAS_FLAT      0
  50. #define BAS_SMOOTH    1
  51.  
  52. /* max number of sides of cylinder to make */
  53. #define BAS_MAXCYLSIDES  512
  54.  
  55. #define PICSTICKBOXWID 185 
  56.  
  57. static ObjPtr AddSticksControls(stick, panelContents)
  58. ObjPtr stick, panelContents;
  59. /*Adds controls to a stick stick*/
  60. {
  61.   ObjPtr var;
  62.   ObjPtr titleBox, button, radio;
  63.   ObjPtr corral, textBox, icon, defaultIcon, name;
  64.   ObjPtr stickField, mainDataset;
  65.   int left, top, bottom, right;
  66.  
  67.   /*Get the stick field*/
  68.   /* The difference between GetObjectVar and GetVar ?? HWU */
  69.   stickField = GetObjectVar("AddSticksControls", stick, MAINDATASET);
  70.   if (!stickField) return ObjFalse;
  71.   while (mainDataset = GetVar(stickField, MAINDATASET))
  72.     {
  73.       stickField = mainDataset;
  74.     }
  75.  
  76.   right = CWINWIDTH - (2 * CORRALBORDER + CWINCORRALWIDTH);
  77.   left = right - PICSTICKBOXWID;
  78.   top = CWINHEIGHT - MINORBORDER;
  79.   bottom = top - (TITLEBOXTOP + 2 * MINORBORDER + 2 * CHECKBOXSPACING + 3 * CHECKBOXHEIGHT);
  80.  
  81.   /* title box for stick style*/
  82.   titleBox = NewTitleBox(left, right, bottom, top,  "Stick Style");
  83.   PrefixList(panelContents, titleBox);
  84.   SetVar(titleBox, PARENT, panelContents);
  85.  
  86.   /* add radio */
  87.   radio = NewRadioButtonGroup("Stick Style Group");
  88.   PrefixList(panelContents, radio);
  89.   SetVar(radio, PARENT, panelContents);
  90.   top -= TITLEBOXTOP + MINORBORDER;
  91.  
  92.   /* Make the buttons */
  93.   button = NewRadioButton(left + MINORBORDER, (left + right) / 2,
  94.               top - CHECKBOXHEIGHT, top, "Lines");
  95.   AddRadioButton(radio, button);
  96.  
  97.   button = NewRadioButton(left + MINORBORDER, (left + right) / 2,
  98.               top - 2 * CHECKBOXHEIGHT - CHECKBOXSPACING, 
  99.               top - CHECKBOXHEIGHT - CHECKBOXSPACING, "Girders");
  100.   AddRadioButton(radio, button);
  101.  
  102.   button = NewRadioButton(left + MINORBORDER, (left + right) / 2,
  103.               top - 3 * CHECKBOXHEIGHT - 2 * CHECKBOXSPACING, 
  104.               top - 2 * CHECKBOXHEIGHT - 2 * CHECKBOXSPACING, 
  105.               "Cylinders");
  106.   AddRadioButton(radio, button);
  107.  
  108.   /* figure out the stickStyle first */
  109.   var = GetIntVar("AddSticksControls", stick, STYLE);
  110.   if (!var)
  111.     {
  112.       SetVar(stick, STYLE, NewInt(BAS_LINES));
  113.     }
  114.  
  115.   AssocDirectControlWithVar(radio, stick, STYLE);
  116.   SetVar(radio, HELPSTRING, NewString("This controls the style of the sticks to be visualized."));
  117.  
  118.   /*Put in the stick corral at the top left*/
  119.   left = MAJORBORDER;
  120.   top = MAJORBORDER;
  121.   corral = NewIconCorral(NULLOBJ, left, left + ONECORRALWIDTH, 
  122.                   CWINHEIGHT - MAJORBORDER - ONECORRALHEIGHT, 
  123.                   CWINHEIGHT - MAJORBORDER, 0);
  124.   SetVar(corral, SINGLECORRAL, ObjTrue);
  125.   SetVar(corral, TOPDOWN, ObjTrue);
  126.   SetVar(corral, NAME, NewString("Stick Field"));
  127.   SetVar(corral, HELPSTRING, NewString("This corral shows the dataset that is \
  128. being used to make the stick display."));
  129.   PrefixList(panelContents, corral);
  130.   SetVar(corral, PARENT, panelContents);
  131.   SetVar(corral, REPOBJ, stick);
  132.   SetMethod(corral, DROPINCONTENTS, DropInMainDatasetCorral);
  133.  
  134.   /*Create the stick source text box*/
  135.   textBox = NewTextBox(left, left + ONECORRALWIDTH, 
  136.      CWINHEIGHT - MAJORBORDER - ONECORRALHEIGHT - TEXTBOXSEP - TEXTBOXHEIGHT,
  137.              CWINHEIGHT - MAJORBORDER - ONECORRALHEIGHT - TEXTBOXSEP,
  138.              0, "Stick Field Text", "Stick Field");
  139.   PrefixList(panelContents, textBox);
  140.   SetVar(textBox, PARENT, panelContents);
  141.   SetTextAlign(textBox, CENTERALIGN);
  142.  
  143.   /*Put in an icon that represents the field*/
  144.   name = GetVar(stickField, NAME);
  145.   defaultIcon = GetVar(stickField, DEFAULTICON);
  146.   if (defaultIcon)
  147.     {
  148.       icon = NewObject(defaultIcon, 0);
  149.       SetVar(icon, NAME, name);
  150.     }
  151.   else
  152.     {
  153.       icon = NewIcon(0, 0, ICONQUESTION, GetString(name));
  154.     }
  155.   SetVar(icon, ICONLOC, NULLOBJ);
  156.   SetVar(icon, REPOBJ, stickField);
  157.   DropIconInCorral(corral, icon);
  158.  
  159.   return ObjTrue;
  160. }
  161.  
  162. #define BAS_I_THRESHOLD 0.8
  163.  
  164. /* chord extension factor */
  165. extern double cefs_2[];
  166.  
  167. static ObjPtr MakeSticksSurface(sticks)
  168. ObjPtr sticks;
  169. /*Makes the surface of a stick object*/
  170. {
  171.     ObjPtr sticksData;
  172.     ObjPtr surface;
  173.     ObjPtr var;
  174.     ObjPtr dataForm;
  175.     ObjPtr edges;            /* edges defining sticks */
  176.     ObjPtr sizeObj;
  177.     int nComponents;         /* # components of data form */
  178.     long nVertices;          /* number of vertices */
  179.     VertexPtr *vertices;     /* all the vertices */
  180.     TwoReals *edgeElements;
  181.     long maxK;               /* number of edges */
  182.     long k;
  183.     TwoReals *residueElements;
  184.     long nResidues;          /* number of residues */
  185.     int stickStyle;          /* style of sticks to make*/
  186.     int colorShading;        /* color shading of sticks to make */
  187.     real size0, size1;       /* if no sizeobj, then both ends of sticks are of size
  188.            size0 ,if is sizeobj, then one end is size0 and the other size1 */
  189.     real sizefactor, sizeoffset, sizeconstant;
  190.     Bool sameForm, capEnds;
  191.  
  192.     /* Get the main dataset*/
  193.     MakeVar(sticks, MAINDATASET);
  194.     if (!(sticksData = GetVar(sticks, MAINDATASET)))
  195.       {
  196.     return ObjFalse;
  197.       }
  198.  
  199.     /* Get its edges */
  200.     if (!(edges = GetDatasetKEdges(sticksData, 1)))
  201.       {
  202.     return ObjFalse;
  203.       }
  204.     /*Get the number of edges */
  205.     maxK = DIMS(edges)[0];
  206.     /* Get edge's tworeal array */
  207.     edgeElements = ELEMENTS(edges);
  208.  
  209.     SetCurForm(FORMFIELD, sticksData);
  210.  
  211.     /* Get number of components, it should be 2 or 3 */
  212.     nComponents = GetNComponents(FORMFIELD);
  213.     if (nComponents != 2 && nComponents != 3)
  214.       {
  215.     return ObjFalse;
  216.       }
  217.  
  218.     /* Get the number of vertices from the current form field */
  219.     if ( 1 != CountTraversalDims(FORMFIELD) )
  220.       {
  221.     return ObjFalse;
  222.       }
  223.     else 
  224.       {
  225.     GetTraversalDims(FORMFIELD, &nVertices);
  226.     if (!nVertices)
  227.       {
  228.         return ObjTrue;
  229.       }
  230.       }
  231.  
  232.     /* Create a new surface picture */
  233.     if (!(surface = NewPicture()))
  234.       {
  235.     return ObjFalse;
  236.       }
  237.  
  238.     /* Caching all the vertices */
  239.     if (!(vertices = (VertexPtr *) Alloc(nVertices * sizeof(VertexPtr))))
  240.       {
  241.     OMErr();
  242.     return ObjFalse;
  243.       }
  244.  
  245.     vertices[0] = NewVertex(surface, VF_FIRSTCANON);
  246.     for (k = 1; k < nVertices; k++)
  247.        {
  248.      vertices[k] = NewVertex(surface, VF_NEXTCANON);
  249.        }
  250.     for (k = 0; k < nVertices; k++)
  251.        {
  252.      vertices[k]->position[0] = SELECTFIELDCOMPONENT(FORMFIELD, 0, &k);
  253.      vertices[k]->position[1] = SELECTFIELDCOMPONENT(FORMFIELD, 1, &k);
  254.      vertices[k]->position[2] = nComponents > 2 ? 
  255.                 SELECTFIELDCOMPONENT(FORMFIELD, 2, &k) : 0.0;
  256.        }
  257.  
  258.     /* Get color shading, 0 for flat color and 1 for smooth color */
  259.     MakeVar(sticks, COLORSHADING);
  260.     var = GetVar(sticks, COLORSHADING);
  261.     if (var)
  262.       {
  263.     colorShading = GetInt(var);
  264.       }
  265.     else
  266.       {
  267.     colorShading = BAS_FLAT; 
  268.       }
  269.  
  270.     /* Get stick style from STYLE variable */
  271.     MakeVar(sticks, STYLE);
  272.     var = GetVar(sticks, STYLE);
  273.     if (var)
  274.       {
  275.     stickStyle = GetInt(var);
  276.       }
  277.     else
  278.       {
  279.     stickStyle = BAS_LINES;
  280.       }
  281.  
  282.     if ((stickStyle == BAS_GIRDERS) || (stickStyle == BAS_CYLINDERS))
  283.       {
  284.         MakeVar(sticks, SIZEFACTOR);
  285.     var = GetRealVar("MakeSticksSurface", sticks, SIZEFACTOR);
  286.     if (var)
  287.       {
  288.         sizefactor = GetReal(var);
  289.       }
  290.         else
  291.       {
  292.         sizefactor = 1.0;
  293.       }
  294.  
  295.         MakeVar(sticks, SIZEOFFSET);
  296.     var = GetRealVar("MakeSticksSurface", sticks, SIZEOFFSET);
  297.     if (var)
  298.       {
  299.         sizeoffset = GetReal(var);
  300.       }
  301.         else
  302.       {
  303.         sizeoffset = 0.0;
  304.       }
  305.  
  306.         sizeObj = NULLOBJ;
  307.     MakeVar(sticks, SIZESWITCH);
  308.     if (GetPredicate(sticks, SIZESWITCH))
  309.       {
  310.             MakeVar(sticks, SIZEOBJ);
  311.         sizeObj = GetVar(sticks, SIZEOBJ);
  312.       }
  313.     if (sizeObj)
  314.       {
  315.          SetCurField(FIELD1, sizeObj);
  316.          SetCurForm(FIELD2, sizeObj);
  317.          sameForm = IdenticalFields(FORMFIELD, FIELD2);
  318.       }
  319.         else
  320.       {
  321.             MakeVar(sticks, SIZECONSTANT);
  322.         var = GetRealVar("MakeSticksSurface", sticks, SIZECONSTANT);
  323.         if (var)
  324.           {
  325.             sizeconstant = GetReal(var);
  326.           }
  327.             else
  328.           {
  329.             sizeconstant = 0.1;
  330.           }
  331.             size0 = sizeconstant * sizefactor + sizeoffset;
  332.       }
  333.  
  334.         /* see if cap is needed by users */
  335.     capEnds=GetPredicate(sticks, CAPENDSP);
  336.       }
  337.  
  338.     if ((stickStyle == BAS_LINES))
  339.       {
  340.     if (colorShading == BAS_SMOOTH)
  341.       {
  342.         VertexPtr line[2];
  343.  
  344.             /* appending all the edges to the surface */
  345.             for (k = 0; k < maxK; k++)
  346.               {
  347.                 line[0] = vertices[(long) edgeElements[k][0]]; 
  348.                 line[1] = vertices[(long) edgeElements[k][1]]; 
  349.         
  350.                 AppendSPolylineToPicture(surface, 1, 0, 2, line);
  351.               }
  352.       }
  353.         else /* colorShading is BAS_FLAT */
  354.       {
  355.         VertexPtr line1[2], line2[2];
  356.  
  357.             for (k = 0; k < maxK; k++)
  358.               {
  359.         
  360.                 line1[0] = vertices[(long) edgeElements[k][0]]; 
  361.                 line2[0] = vertices[(long) edgeElements[k][1]]; 
  362.  
  363.         line1[1] = InsertVertex(surface, line1[0]);
  364.         line2[1] = InsertVertex(surface, line2[0]);
  365.  
  366.         line1[1]->position[0] = line2[1]->position[0] = 
  367.            0.5 * (line1[0]->position[0] + line2[0]->position[0]);
  368.         line1[1]->position[1] = line2[1]->position[1] = 
  369.            0.5 * (line1[0]->position[1] + line2[0]->position[1]);
  370.         line1[1]->position[2] = line2[1]->position[2] = 
  371.            0.5 * (line1[0]->position[2] + line2[0]->position[2]);
  372.  
  373.                 AppendSPolylineToPicture(surface, 1, 0, 2, line1);
  374.                 AppendSPolylineToPicture(surface, 1, 0, 2, line2);
  375.               }
  376.       }
  377.       }
  378.     else if ((stickStyle == BAS_GIRDERS))
  379.       {
  380.     VertexPtr edgeVertices[2];
  381.     float a[3];          /* axial vectors */
  382.     float r[4][3];       /* the normalized radial vectors */
  383.     float rSized0[4][3], rSized1[4][3];  /* sized radial vectors */
  384.     VertexPtr gutterEnd0[8], gutterEnd1[8];
  385.     PolysPtr polys;
  386.     int i, m;
  387.  
  388.     polys = AppendPolysToPicture(surface);
  389.  
  390.         for (k = 0; k < maxK; k++)
  391.            {
  392.              edgeVertices[0] = vertices[(long) edgeElements[k][0]]; 
  393.              edgeVertices[1] = vertices[(long) edgeElements[k][1]]; 
  394.   
  395.              /* axial vector */
  396.          for(i = 0; i < 3; i++)
  397.            {
  398.              a[i] = edgeVertices[0]->position[i] - edgeVertices[1]->position[i];
  399.            }
  400.          NORMALIZE(a);
  401.   
  402.          /* get first radial vector */
  403.              if (ABS(a[0]) >= BAS_I_THRESHOLD)
  404.            {
  405.              /* by crossing with j axis vector (0, 1, 0) */
  406.              r[0][0] = -a[2]; r[0][1] = 0.0; r[0][2] = a[0]; 
  407.            }
  408.              else
  409.            {
  410.              /* by crossing with i axis vector (1, 0, 0) */
  411.              r[0][0] = 0.0; r[0][1] = a[2]; r[0][2] = -a[1]; 
  412.            }
  413.              /* needs to normalize r[0] */
  414.          NORMALIZE(r[0]);
  415.             
  416.          /* get orgothonal radial vector, counterclockwisely */
  417.          CROSS(a, r[0], r[1]);
  418.          /* do not need to normalize r[1] since r[0] and a are both unit vector and they are perpendicular to each other, r[1] must be unit vector also*/
  419.  
  420.              /* consider the size factor and calculate their opposite direction*/
  421.          /* r[0] and r[2] are opposite, and r[1] and r[3] opposite*/
  422.          if (sizeObj)
  423.            {
  424.          long index0, index1;
  425.          real sizeval0, sizeval1;
  426.  
  427.          index0 = (long) edgeElements[k][0];
  428.          index1 = (long) edgeElements[k][1];
  429.          if (sameForm)
  430.            {
  431.              sizeval0 = SelectFieldScalar(FIELD1, &index0);
  432.              sizeval1 = SelectFieldScalar(FIELD1, &index1);
  433.            }
  434.                  else
  435.            {
  436.              sizeval0 = SampleSpatScalar(FIELD1, FIELD2, 3, 
  437.                          edgeVertices[0]->position, true);
  438.              sizeval1 = SampleSpatScalar(FIELD1, FIELD2, 3, 
  439.                          edgeVertices[1]->position, true);
  440.            }
  441.          if (sizeval0 == missingData)
  442.            {
  443.              sizeval0 = 1.0;
  444.            }
  445.          if (sizeval1 == missingData)
  446.            {
  447.              sizeval1 = 1.0;
  448.            }
  449.          size0 = sizeval0 * sizefactor + sizeoffset;
  450.          size1 = sizeval1 * sizefactor + sizeoffset;
  451.  
  452.              for (m = 0; m < 3; m++)
  453.                 {
  454.                   r[2][m] = -r[0][m];
  455.               rSized0[0][m] = r[0][m] * size0;
  456.               rSized1[0][m] = r[0][m] * size1;
  457.               rSized0[2][m] = -rSized0[0][m];
  458.               rSized1[2][m] = -rSized1[0][m];
  459.               r[3][m] = -r[1][m];
  460.               rSized0[1][m] = r[1][m] * size0;
  461.               rSized1[1][m] = r[1][m] * size1;
  462.               rSized0[3][m] = -rSized0[1][m];
  463.               rSized1[3][m] = -rSized1[1][m];
  464.             } 
  465.            }
  466.              else /* no sizeObj */
  467.            {
  468.              for (m = 0; m < 3; m++)
  469.                 {
  470.                   r[2][m] = -r[0][m];
  471.               rSized0[0][m] = r[0][m] * size0;
  472.               rSized0[2][m] = -rSized0[0][m];
  473.               r[3][m] = -r[1][m];
  474.               rSized0[1][m] = r[1][m] * size0;
  475.               rSized0[3][m] = -rSized0[1][m];
  476.             } 
  477.            }
  478.   
  479.          /* vertex position */
  480.          /* group as (0,1) (2,3) (4,5) (6,7) as of the same position */
  481.          for (i = 0; i < 4; i++)
  482.             {
  483.                   gutterEnd0[2 * i] =     InsertVertex(surface, edgeVertices[0]);
  484.                   gutterEnd0[2 * i + 1] = InsertVertex(surface, edgeVertices[0]);
  485.                   gutterEnd1[2 * i] =     InsertVertex(surface, edgeVertices[1]);
  486.                   gutterEnd1[2 * i + 1] = InsertVertex(surface, edgeVertices[1]);
  487.  
  488.           if (sizeObj)
  489.             {
  490.                   for (m = 0; m < 3; m++)
  491.                      {
  492.                    gutterEnd0[2 * i]->position[m] =
  493.                    gutterEnd0[2 * i + 1]->position[m] =
  494.                      edgeVertices[0]->position[m] + rSized0[i][m];
  495.                    gutterEnd1[2 * i]->position[m] =
  496.                    gutterEnd1[2 * i + 1]->position[m] =
  497.                      edgeVertices[1]->position[m] + rSized1[i][m];
  498.                          }
  499.             }
  500.                   else  /* no sizeObj */
  501.             {
  502.                   for (m = 0; m < 3; m++)
  503.                      {
  504.                    gutterEnd0[2 * i]->position[m] =
  505.                    gutterEnd0[2 * i + 1]->position[m] =
  506.                      edgeVertices[0]->position[m] + rSized0[i][m];
  507.                    gutterEnd1[2 * i]->position[m] =
  508.                    gutterEnd1[2 * i + 1]->position[m] =
  509.                      edgeVertices[1]->position[m] + rSized0[i][m];
  510.                          }
  511.                     }
  512.                 }
  513.              /* vertex normal*/
  514.          if (sizeObj && size0 != size1)
  515.            {
  516.                  /* get normal of the first surface, with vertices 7 and 0 */ 
  517.          float axialv[3], verticalv[3], nv[3], nextnv[3];
  518.  
  519.          for (m = 0; m < 3; m++)
  520.             {
  521.               axialv[m] = gutterEnd0[0]->position[m] -
  522.                            gutterEnd1[0]->position[m];
  523.                       verticalv[m] = gutterEnd1[7]->position[m] -
  524.                          gutterEnd1[0]->position[m];
  525.             }
  526.                  CROSS(axialv, verticalv, nv);
  527.          NORMALIZE(nv);
  528.          for (m = 0; m < 3; m++)
  529.             {
  530.                       gutterEnd0[7]->normal[m] = 
  531.                       gutterEnd0[0]->normal[m] = 
  532.               gutterEnd1[7]->normal[m] = 
  533.               gutterEnd1[0]->normal[m] = nv[m]; 
  534.             }
  535.  
  536.                  /* now get the normals for the other surface by crossing
  537.          its previous surface's normal with the boundary vector*/
  538.          for (i = 1; i < 4; i++)
  539.             {
  540.               for (m = 0; m < 3; m++)
  541.                  {
  542.                axialv[m] = gutterEnd0[2 * i - 1]->position[m] -
  543.                            gutterEnd1[2 * i - 1]->position[m];
  544.                  }
  545.                       CROSS(axialv, nv, nextnv);
  546.               NORMALIZE(nextnv);
  547.               for (m = 0; m < 3; m++)
  548.              {
  549.                            gutterEnd0[2 * i - 1]->normal[m] = 
  550.                            gutterEnd0[2 * i]->normal[m] = 
  551.                    gutterEnd1[2 * i - 1]->normal[m] = 
  552.                    gutterEnd1[2 * i]->normal[m] = nv[m] = nextnv[m]; 
  553.              }
  554.                     }
  555.  
  556.            }
  557.              else /* no sizeObj, or size0 == size 1 */
  558.            {
  559.              for (i = 1; i < 4; i++)
  560.                 for (m = 0; m < 3; m++)
  561.                    {
  562.                          gutterEnd0[2 * i - 1]->normal[m] = 
  563.                          gutterEnd0[2 * i]->normal[m] = 
  564.                  gutterEnd1[2 * i - 1]->normal[m] = 
  565.                  gutterEnd1[2 * i]->normal[m] = 
  566.                  SQ22 * (r[i - 1][m] + r[i][m]);
  567.                }
  568.              for (m = 0; m < 3; m++)
  569.                 {
  570.                       gutterEnd0[7]->normal[m] = 
  571.                       gutterEnd0[0]->normal[m] = 
  572.                   gutterEnd1[7]->normal[m] = 
  573.                   gutterEnd1[0]->normal[m] = SQ22 * (r[3][m] + r[0][m]);
  574.                 }
  575.                }
  576.  
  577.              /* constructing the surface with the vertices */
  578.          if (colorShading == BAS_SMOOTH)
  579.            {
  580.              VertexPtr v[4];
  581.  
  582.          for (i = 1; i < 4; i++)
  583.             {
  584.                   v[0] = gutterEnd0[2 * i - 1];
  585.                   v[1] = gutterEnd0[2 * i];
  586.                   v[2] = gutterEnd1[2 * i];
  587.                   v[3] = gutterEnd1[2 * i - 1];
  588.                   AppendSPolyToPolys(polys, 4, v);
  589.             }
  590.                  v[0] = gutterEnd0[7];
  591.                  v[1] = gutterEnd0[0];
  592.                  v[2] = gutterEnd1[0];
  593.                  v[3] = gutterEnd1[7];
  594.                  AppendSPolyToPolys(polys, 4, v);
  595.                }
  596.              else /* colorShading is BAS_FLAT */
  597.            {
  598.              VertexPtr gutterMid0[8], gutterMid1[8];
  599.              VertexPtr v[4];
  600.      
  601.              /* vertex position */
  602.              /* group as (0,1) (2,3) (4,5) (6,7) as of the same position */
  603.              for (i = 0; i < 4; i++)
  604.                 {
  605.                       gutterMid0[2 * i] = InsertVertex(surface, edgeVertices[0]);
  606.                       gutterMid0[2 * i + 1] = InsertVertex(surface, edgeVertices[0]);
  607.                       gutterMid1[2 * i] = InsertVertex(surface, edgeVertices[1]);
  608.                       gutterMid1[2 * i + 1] = InsertVertex(surface, edgeVertices[1]);
  609.  
  610.                   for (m = 0; m < 3; m++)
  611.                      {
  612.                    gutterMid0[2 * i]->position[m] =
  613.                    gutterMid0[2 * i + 1]->position[m] =
  614.                    gutterMid1[2 * i]->position[m] =
  615.                    gutterMid1[2 * i + 1]->position[m] =
  616.                            0.5 * (gutterEnd0[2 * i]->position[m] + 
  617.                        gutterEnd1[2 * i]->position[m]);
  618.                          }
  619.                     }
  620.                  /* vertex normal*/
  621.              for (i = 0; i < 8; i++)
  622.                 for (m = 0; m < 3; m++)
  623.                    {
  624.                          gutterMid0[i]->normal[m] = 
  625.                  gutterMid1[i]->normal[m] = 
  626.                  gutterEnd0[i]->normal[m];
  627.                }
  628.                  for (i = 1; i < 4; i++)
  629.             {
  630.                       v[0] = gutterEnd0[2 * i - 1];
  631.                   v[1] = gutterEnd0[2 * i];
  632.                   v[2] = gutterMid0[2 * i];
  633.                   v[3] = gutterMid0[2 * i - 1];
  634.                   AppendSPolyToPolys(polys, 4, v);
  635.      
  636.                       v[0] = gutterEnd1[2 * i - 1];
  637.                       v[1] = gutterEnd1[2 * i];
  638.                       v[2] = gutterMid1[2 * i];
  639.                       v[3] = gutterMid1[2 * i - 1];
  640.                       AppendSPolyToPolys(polys, 4, v);
  641.             }
  642.                  v[0] = gutterEnd0[7];
  643.                  v[1] = gutterEnd0[0];
  644.                  v[2] = gutterMid0[0];
  645.                  v[3] = gutterMid0[7];
  646.                  AppendSPolyToPolys(polys, 4, v);
  647.      
  648.                  v[0] = gutterEnd1[7];
  649.                  v[1] = gutterEnd1[0];
  650.                  v[2] = gutterMid1[0];
  651.                  v[3] = gutterMid1[7];
  652.                  AppendSPolyToPolys(polys, 4, v);
  653.            } /* else */
  654.  
  655.              /* end caps */
  656.          if (capEnds)
  657.            {
  658.                  VertexPtr cap0[4], cap1[4];
  659.  
  660.              for (i = 0; i < 4; i++)
  661.                 {
  662.                       cap0[i] = InsertVertex(surface, edgeVertices[0]);
  663.                       cap1[i] = InsertVertex(surface, edgeVertices[1]);
  664.  
  665.               for (m = 0; m < 3; m++)
  666.                  {
  667.                    cap0[i]->position[m]=gutterEnd0[2 * i]->position[m];
  668.                    cap1[i]->position[m]=gutterEnd1[2 * i]->position[m];
  669.                            cap0[i]->normal[m] = a[m];
  670.                            cap1[i]->normal[m] = - a[m];
  671.                          }
  672.                 }
  673.              AppendSPolyToPolys(polys, 4, cap0);
  674.              AppendSPolyToPolys(polys, 4, cap1);
  675.            }
  676.           } /* for loop */
  677.       }
  678.     else if ((stickStyle == BAS_CYLINDERS))
  679.       {
  680.     VertexPtr edgeVertices[2];
  681.     VertexPtr cylEnd0[BAS_MAXCYLSIDES], cylEnd1[BAS_MAXCYLSIDES]; 
  682.     float a[3];          /* axial vectors */
  683.     float r[BAS_MAXCYLSIDES][3];       /* the normalized radial vectors */
  684.     float rSized0[BAS_MAXCYLSIDES][3], 
  685.           rSized1[BAS_MAXCYLSIDES][3]; /* if size are diff.  */
  686.     PolysPtr polys;
  687.     int nSubdivisions;
  688.     int curNSides = 4;
  689.     int delta;
  690.     int i, m;
  691.  
  692.         /* Get number of subdivisions the users want */
  693.         MakeVar(sticks, FRUSTUMSUBDIV);
  694.         var = GetVar(sticks, FRUSTUMSUBDIV);
  695.         if (var)
  696.           {
  697.         nSubdivisions = GetInt(var);
  698.           }
  699.         else
  700.           {
  701.         nSubdivisions = 2;  /* Octal sides */
  702.           }
  703.         while (nSubdivisions--) curNSides *= 2;
  704.  
  705.     polys = AppendPolysToPicture(surface);
  706.  
  707.         for (k = 0; k < maxK; k++)
  708.            {
  709.              edgeVertices[0] = vertices[(long) edgeElements[k][0]]; 
  710.              edgeVertices[1] = vertices[(long) edgeElements[k][1]]; 
  711.   
  712.          for (i = 0; i < curNSides; i++)
  713.             {
  714.                   cylEnd0[i] = InsertVertex(surface, edgeVertices[0]);
  715.                   cylEnd1[i] = InsertVertex(surface, edgeVertices[1]);
  716.             }
  717.              /* axial vector */
  718.          for (i = 0; i < 3; i++)
  719.             {
  720.               a[i] = edgeVertices[0]->position[i] - edgeVertices[1]->position[i];
  721.             }
  722.              /* a is the normal vector of the plane at e1  */
  723.          NORMALIZE(a);
  724.   
  725.          /* get first radial vector */
  726.              if (ABS(a[0]) >= BAS_I_THRESHOLD)
  727.            {
  728.              /* by crossing with j axis vector (0, 1, 0) */
  729.              r[0][0] = -a[2]; r[0][1] = 0.0; r[0][2] = a[0]; 
  730.            }
  731.              else
  732.            {
  733.              /* by crossing with i axis vector (1, 0, 0) */
  734.              r[0][0] = 0.0; r[0][1] = a[2]; r[0][2] = -a[1]; 
  735.            }
  736.              /* needs to normalize r[0] */
  737.          NORMALIZE(r[0]);
  738.             
  739.          delta = curNSides / 4;
  740.          /* get 3 other main radial vectors */
  741.          CROSS(a, r[0], r[delta]);
  742.              /* consider the size factor and calculate their opposite direction*/
  743.          /* r[0] and r[2] are opposite, and r[1] and r[3] opposite*/
  744.          if (sizeObj)
  745.            {
  746.          long index0, index1;
  747.          real sizeval0, sizeval1;
  748.  
  749.          index0 = (long) edgeElements[k][0];
  750.          index1 = (long) edgeElements[k][1];
  751.          if (sameForm)
  752.            {
  753.              sizeval0 = SelectFieldScalar(FIELD1, &index0);
  754.              sizeval1 = SelectFieldScalar(FIELD1, &index1);
  755.            }
  756.                  else
  757.            {
  758.              sizeval0 = SampleSpatScalar(FIELD1, FIELD2, 3, 
  759.                          edgeVertices[0]->position, true);
  760.              sizeval1 = SampleSpatScalar(FIELD1, FIELD2, 3, 
  761.                          edgeVertices[1]->position, true);
  762.            }
  763.          if (sizeval0 == missingData)
  764.            {
  765.              sizeval0 = 1.0;
  766.            }
  767.          if (sizeval1 == missingData)
  768.            {
  769.              sizeval1 = 1.0;
  770.            }
  771.          size0 = sizeval0 * sizefactor + sizeoffset;
  772.          size1 = sizeval1 * sizefactor + sizeoffset;
  773.  
  774.              for (m = 0; m < 3; m++)
  775.                 {
  776.                   r[2 * delta][m] = -r[0][m];
  777.               rSized0[0][m] = r[0][m] * size0;
  778.               rSized1[0][m] = r[0][m] * size1;
  779.               rSized0[2 * delta][m] = -rSized0[0][m];
  780.               rSized1[2 * delta][m] = -rSized1[0][m];
  781.               r[3 * delta][m] = -r[delta][m];
  782.               rSized0[delta][m] = r[delta][m] * size0;
  783.               rSized1[delta][m] = r[delta][m] * size1;
  784.               rSized0[3 * delta][m] = -rSized0[delta][m];
  785.               rSized1[3 * delta][m] = -rSized1[delta][m];
  786.             } 
  787.            }
  788.              else /* no sizeObj */
  789.            {
  790.              for (m = 0; m < 3; m++)
  791.                 {
  792.                   r[2*delta][m] = -r[0][m];
  793.               rSized0[0][m] = r[0][m] * size0;
  794.               rSized0[2 * delta][m] = -rSized0[0][m];
  795.               r[3 * delta][m] = -r[delta][m];
  796.               rSized0[delta][m] = r[delta][m] * size0;
  797.               rSized0[3 * delta][m] = -rSized0[delta][m];
  798.             } 
  799.            }
  800.              /* get the remaining radial vectors */
  801.          i = 0;
  802.          for ( ;delta > 1; delta /= 2 )
  803.         {
  804.           register double cef_2;
  805.           int beg;
  806.  
  807.           cef_2 = cefs_2[i];
  808.           i++;
  809.  
  810.           for (beg = 0; beg < curNSides; beg += delta)
  811.              {
  812.                register int mid, end; /* mid and end points */
  813.  
  814.                end = (beg + delta) & (curNSides - 1); /* circular */
  815.                mid = beg + delta / 2;
  816.  
  817.                for(m = 0; m < 3; m++)
  818.              {
  819.                r[mid][m] = cef_2 * (r[beg][m] + r[end][m]);
  820.                rSized0[mid][m] = r[mid][m] * size0;
  821.                if (sizeObj)
  822.                  {
  823.                    rSized1[mid][m] = r[mid][m] * size1;
  824.                  }
  825.              }
  826.              }
  827.         }
  828.   
  829.          delta = curNSides / 4;
  830.          for (i = 0; i < curNSides; i++)
  831.         {
  832.               for (m = 0; m < 3; m++)
  833.                  {
  834.                    /* vertex position */
  835.                cylEnd0[i]->position[m] =
  836.                  edgeVertices[0]->position[m] + rSized0[i][m];
  837.                cylEnd1[i]->position[m] = (sizeObj ? 
  838.                  edgeVertices[1]->position[m] + rSized1[i][m] :
  839.                  edgeVertices[1]->position[m] + rSized0[i][m]);
  840.                      }
  841.             
  842.                   /* vetex normal */
  843.           if (sizeObj && size0 != size1)
  844.             {
  845.               int j;
  846.               real boundaryVector[3], nv[3];
  847.   
  848.               /* circular clipping */
  849.               j = ( i + delta ) & (curNSides -1);
  850.               for (m = 0; m < 3; m++)
  851.              {
  852.                boundaryVector[m] = cylEnd0[i]->position[m] -
  853.                            cylEnd1[i]->position[m];
  854.              }
  855.                       CROSS(r[j], boundaryVector, nv);
  856.               NORMALIZE(nv);
  857.               for (m = 0; m < 3; m++)
  858.                  {
  859.                            cylEnd0[i]->normal[m] = cylEnd1[i]->normal[m] = nv[m];
  860.                          }
  861.             }
  862.                   else
  863.             {
  864.               for (m = 0; m < 3; m++)
  865.              {
  866.                            cylEnd0[i]->normal[m] = cylEnd1[i]->normal[m] = r[i][m];
  867.                          }
  868.                     }
  869.                  }
  870.  
  871.          if (colorShading == BAS_SMOOTH)
  872.            {
  873.              VertexPtr v[4];
  874.  
  875.                  /* for cylinder body */
  876.          for (i = 0; i < curNSides - 1; i++)
  877.             {
  878.                   v[0] = cylEnd0[i];
  879.                   v[1] = cylEnd0[i + 1];
  880.                   v[2] = cylEnd1[i + 1];
  881.                   v[3] = cylEnd1[i];
  882.                   AppendSPolyToPolys(polys, 4, v);
  883.             }
  884.              v[0] = cylEnd0[curNSides - 1];
  885.              v[1] = cylEnd0[0];
  886.              v[2] = cylEnd1[0];
  887.              v[3] = cylEnd1[curNSides - 1];
  888.              AppendSPolyToPolys(polys, 4, v);
  889.                }
  890.              else  /* colorShading is BAS_FLAT */
  891.            {
  892.              VertexPtr cylMid0[BAS_MAXCYLSIDES], cylMid1[BAS_MAXCYLSIDES]; 
  893.          VertexPtr v[4];
  894.  
  895.              for (i = 0; i < curNSides; i++)
  896.             {
  897.                       cylMid0[i] = InsertVertex(surface, edgeVertices[0]);
  898.                       cylMid1[i] = InsertVertex(surface, edgeVertices[1]);
  899.             }
  900.              for (i = 0; i < curNSides; i++)
  901.                 for (m = 0; m < 3; m++)
  902.                {
  903.                      /* vertex position */
  904.                  cylMid0[i]->position[m] = cylMid1[i]->position[m] =
  905.              0.5 * (cylEnd0[i]->position[m] + cylEnd1[i]->position[m]);
  906.                          /* vertex normal */
  907.                          cylMid0[i]->normal[m] = cylMid1[i]->normal[m] 
  908.                                    = cylEnd0[i]->normal[m];
  909.                        }
  910.  
  911.                  /* constructing the surface with the vertices */
  912.          for (i = 0; i < curNSides - 1; i++)
  913.             {
  914.                   v[0] = cylEnd0[i];
  915.                   v[1] = cylEnd0[i + 1];
  916.                   v[2] = cylMid0[i + 1];
  917.                   v[3] = cylMid0[i];
  918.                   AppendSPolyToPolys(polys, 4, v);
  919.  
  920.                   v[0] = cylEnd1[i];
  921.                   v[1] = cylEnd1[i + 1];
  922.                   v[2] = cylMid1[i + 1];
  923.                   v[3] = cylMid1[i];
  924.                   AppendSPolyToPolys(polys, 4, v);
  925.            }
  926.              v[0] = cylEnd0[curNSides - 1];
  927.              v[1] = cylEnd0[0];
  928.              v[2] = cylMid0[0];
  929.              v[3] = cylMid0[curNSides - 1];
  930.              AppendSPolyToPolys(polys, 4, v);
  931.  
  932.              v[0] = cylEnd1[curNSides - 1];
  933.              v[1] = cylEnd1[0];
  934.              v[2] = cylMid1[0];
  935.              v[3] = cylMid1[curNSides - 1];
  936.              AppendSPolyToPolys(polys, 4, v);
  937.            }
  938.  
  939.              /* for end caps */
  940.               if (capEnds)
  941.         {
  942.                   VertexPtr cap0[BAS_MAXCYLSIDES], cap1[BAS_MAXCYLSIDES];
  943.  
  944.               for (i = 0; i < curNSides; i++)
  945.              {
  946.                        cap0[i] = InsertVertex(surface, edgeVertices[0]);
  947.                        cap1[i] = InsertVertex(surface, edgeVertices[1]);
  948.                    for (m = 0; m < 3; m++)
  949.                   {
  950.                 cap0[i]->position[m] = cylEnd0[i]->position[m];
  951.                 cap1[i]->position[m] = cylEnd1[i]->position[m];
  952.                             cap0[i]->normal[m] = a[m];
  953.                             cap1[i]->normal[m] = - a[m];
  954.                   }
  955.              }
  956.  
  957.               AppendSPolyToPolys(polys, curNSides, cap0);
  958.               AppendSPolyToPolys(polys, curNSides, cap1);
  959.                 }
  960.         } /* for loop */
  961.       }
  962.     else
  963.       {
  964.     /* no surface style */
  965.       }
  966.  
  967.     /* set the pictures's MAINDATASET for facilitating ColorPictureByObject */
  968.     SetVar(surface, MAINDATASET, sticksData);
  969.  
  970.     SetVar(sticks, SURFACE, surface);
  971.     SetVar(surface, REPOBJ, sticks);
  972.     Free(vertices);
  973.     return ObjTrue;
  974. }
  975.  
  976. static ObjPtr SticksInit(sticks)
  977. ObjPtr sticks;
  978. /*Initializes a stick object*/
  979. {
  980.     ObjPtr colorObj;
  981.  
  982.     /* here we make maindataset and set up colorObj for sticks */
  983.     MakeVar(sticks, MAINDATASET);
  984.     SetVar(sticks, COLOROBJ, colorObj = GetVar(sticks, MAINDATASET));
  985.  
  986.     if (colorObj)
  987.     {
  988.       ObjPtr mainDataset, sizeObj;
  989.  
  990.       /* set the stick vis object as having color attribute */
  991.       SetVar(sticks, COLORS, ObjTrue);
  992.  
  993.       /* see if there's a sizeobj associated with its main dataset */
  994.       while (mainDataset = GetVar(colorObj, MAINDATASET))
  995.      {
  996.        colorObj = mainDataset;
  997.      }
  998.       sizeObj = GetVar(colorObj, SIZEOBJ);
  999.       if (sizeObj)
  1000.     {
  1001.       SetVar(sticks, SIZEOBJ, sizeObj);
  1002.       SetVar(sticks, SIZESWITCH, NewInt(1));
  1003.     }
  1004.     }
  1005.  
  1006.     return ObjTrue;
  1007. }
  1008.  
  1009. static ObjPtr SetSticksMainDataset(visObj, dataSet)
  1010. ObjPtr visObj, dataSet;
  1011. /*Sets the main data set of visObj to dataSet*/
  1012. {
  1013.     SetVar(visObj, MAINDATASET, dataSet);
  1014.     return ObjTrue;
  1015. }
  1016.  
  1017. /* Initializes sticks visualization object */
  1018. void InitSticks()
  1019. {
  1020.     ObjPtr icon;
  1021.  
  1022.     /*Class for sticks */
  1023.     visSticks = NewObject(visSized, 0);
  1024.     AddToReferenceList(visSticks);
  1025.  
  1026.     SetVar(visSticks, NAME, NewString("Sticks"));
  1027.     SetVar(visSticks, STYLE, NewInt(BAS_LINES)); 
  1028.     SetVar(visSticks, SIZECONSTANT, NewReal(0.1));
  1029.     SetVar(visSticks, SIZEFACTOR, NewReal(1.0));
  1030.     SetVar(visSticks, SIZEOFFSET, NewReal(0.0));
  1031.  
  1032.     SetMethod(visSticks, INITIALIZE, SticksInit);
  1033.  
  1034.     /* notice a new icon is created for this stick visobj */
  1035.     SetVar(visSticks, DEFAULTICON, icon = NewObject(visIcon, 0));
  1036.     SetVar(icon, WHICHICON, NewInt(ICONSTICKS));
  1037.     SetVar(icon, NAME, NewString("Sticks"));
  1038.     SetVar(icon, HELPSTRING,
  1039.     NewString("This icon represents a stick visualization object.  \
  1040. This object is used to display a skeleton of molecular structure."));
  1041.  
  1042.     DeclareIndirectDependency(visSticks, SURFACE, MAINDATASET, DATA);
  1043.     DeclareDependency(visSticks, SURFACE, STYLE);
  1044.     /* surface depends on colorshading */
  1045.     SetVar(visSticks, COLORSHADING, ObjFalse); /* set colorshading flat */
  1046.     DeclareDependency(visSticks, SURFACE, COLORSHADING);
  1047.  
  1048.     SetMethod(visSticks, SURFACE, MakeSticksSurface);
  1049.     SetMethod(visSticks, SETMAINDATASET, SetSticksMainDataset);
  1050.     SetMethod(visSticks, ADDCONTROLS, AddSticksControls);
  1051.  
  1052.     icon = NewIcon(0, 0, ICONSTICKS, "Sticks");
  1053.     SetVar(icon, HELPSTRING,
  1054.     NewString("Click on this icon to see a panel of controls for the stick object."));
  1055.     SetVar(visSticks, CONTROLICON, icon);
  1056.  
  1057.     DefineVisMapping(DS_HASFORM | DS_HASFIELD | DS_UNSTRUCTURED, 1, 3, 1, visSticks);
  1058.     DefineVisMapping(DS_HASFORM | DS_HASFIELD | DS_UNSTRUCTURED, 1, 2, 1, visSticks);
  1059. }
  1060.  
  1061. void KillSticks()
  1062. /*Kills sticks visualization*/
  1063. {
  1064.     DeleteThing(visSticks);
  1065. }
  1066.